#ifndef _RBFFF_H_
#define _RBFFF_H_

#include "PointSet.h"
#include "OctTree.h"
//#include "math.h"
#include "ImplicitFunction.h"

class RBF : public ImplicitFunction
{
public:
	RBF();
	virtual ~RBF();

	void setPointSet(PointSet *ps, bool takeOwnership);
  void setSubRBF(RBF *rbf);
  void setGlobal() { isGlobal = true; };
  void setSupportConstraint(int minPoints);
  void getDimension(float &x, float &y, float &z);
  //bool getNormal(float x, float y, float z, float delta, float &nx, float &ny, float &nz);

	//Octree oprtations
  const PointSet* getPointSet() { return ps; };
  void splitOctreeChild(int M);
  void computeOctreeCenter();
	float getOctreeAveragedLeafSize();
  float getOctreePeak();
  void generateOctreePointList(int &size, float x, float y, float z, float r, int max) { tree->generatePointList(size, x, y, z, r, max); };
  float getOctreePointX(int index) { return tree->getPointX(index); };
  float getOctreePointY(int index) { return tree->getPointY(index); };
  float getOctreePointZ(int index) { return tree->getPointZ(index); };
  float getOctreeNormalX(int index) { return tree->getNormalX(index); };
  float getOctreeNormalY(int index) { return tree->getNormalY(index); };
  float getOctreeNormalZ(int index) { return tree->getNormalZ(index); };

	float valueAndGradient(float g[3], float x, float y, float z);

	//Fitting
	void fit(float support, float ani, float smooth);
	void computeQuad(float support);

	//Wendland's CSRBF
	void setSupport(float T);

private:
	//Fitting
	void computeSolution(float T);

  //ImplicitFunction virtual functions
	float value(float x, float y, float z);
	int numSupportPoints(float x, float y, float z);

  // Test for boundary
  BoundaryTest testBounds(float x, float y, float z);

	//Wendland's CSRBF
	inline void weightD(float g[], double d2, float vx, float vy, float vz);
	inline double weight(double d2);

	inline void crossProduct(float n[3], const float v1[3], const float v2[3]);
	inline double dotProduct(const float v1[3], const float v2[3]);

private:
	OctTree* tree;
	PointSet* ps;
	RBF* sub_rbf; //Previous level function

  float bounds[6];
  int minSupport;
	double *sol; 	//Shift amount
	float support;
	float c;
	double T2; 	//T2 = support^2
	double K2;
	bool isGlobal; 	//Base function or not?
  bool ownPointSet;
};

#endif
